
Image(filename = 'picture/1.0.1.jpg', width=800, height=400)


# Correlation matrix
Image(filename = 'picture/2.1.1.jpg', width=500, height=500)


# Summary
Image(filename = 'picture/2.1.2.jpg', width=600, height=400)


# Choropleth map
Image(filename = 'picture/2.2.1.jpg', width=800, height=400)


# barlot
Image(filename = 'picture/2.3.1.jpg', width=800, height=400)


# Regression model analysis: Linear vs. Random Forest
Image(filename = 'picture/2.4.1.jpg', width=700, height=300)


# Actual vs Predicted Happiness Scores using Random Forest
Image(filename = 'picture/2.4.2.jpg', width=800, height=400)


# Actual vs Predicted Happiness Scores using Random Forest
Image(filename = 'picture/2.5.1.jpg', width=600, height=400)


# Actual vs Predicted Happiness Scores using Random Forest
Image(filename = 'picture/2.5.2.jpg', width=1000, height=400)


# 1.1 Code
# Import Packages
import numpy as np
import pandas as pd
import altair as alt
import seaborn as sns
import plotly.express as px
import statsmodels.api as sm
import matplotlib.pyplot as plt

from sklearn import preprocessing
from IPython.display import Image # Insert pictures
from sklearn.cluster import KMeans
from sklearn.preprocessing import StandardScaler
from scipy.spatial.distance import cdist #calculate distance

# disable row limit for plotting
alt.data_transformers.disable_max_rows()

# View the WHR Dataset
whr = pd.read_csv('data/whr-2023.csv')
whr.head()


# Code of 2.1
# Selecting relevant columns for the correlation analysis
columns_of_interest = ['Life Ladder', 'Log GDP per capita', 'Social support', 
                       'Healthy life expectancy at birth', 
                       'Freedom to make life choices', 'Positive affect','Generosity', 
                       'Negative affect','Perceptions of corruption']
whr_subset = whr[columns_of_interest]

# Calculating the correlation matrix
corr_matrix = whr_subset.corr()

# Plotting the correlation matrix
plt.figure(figsize=(8, 6))
sns.heatmap(corr_matrix, annot=True, cmap='coolwarm')
plt.title('Correlation Matrix of Happiness Score and Socio-Economic Factors')
plt.show()


# Selecting independent variables
columns_of_interest2 = ['Log GDP per capita', 'Social support', 'Healthy life expectancy at birth', 
                        'Freedom to make life choices', 'Generosity', 'Perceptions of corruption',
                        'Positive affect','Negative affect']

# Handling missing and infinite values
whr_na = whr.replace([np.inf, -np.inf], np.nan)  
whr_na = whr_na.dropna(subset=columns_of_interest2 + ['Life Ladder']) 

X = whr_na[columns_of_interest2]
X = sm.add_constant(X)  # Adding a constant to the model

Y = whr_na['Life Ladder']

# Fitting the model
model = sm.OLS(Y, X).fit()

# Printing the summary
print(model.summary())


# Code of 2.2
# Load the iso datasets
country_codes = pd.read_csv('data/iso-country-codes.csv')

# Find the 'Life Ladder' for the 1st and last available year
whr_first = whr.groupby('Country name').first().reset_index()
whr_last = whr.groupby('Country name').last().reset_index()

# Calculate the difference
whr_diff = pd.DataFrame()
whr_diff['Country name'] = whr_last['Country name']
whr_diff['Ladder Difference'] = whr_last['Life Ladder'] - whr_first['Life Ladder']

# Merge the country codes
merged_data = whr_diff.merge(country_codes, how='left', left_on='Country name', 
                             right_on='English short name lower case')

# Create the choropleth map
choropleth_map = px.choropleth(merged_data,
                    locations='Alpha-3 code',  
                    color='Ladder Difference',  
                    hover_name='Country name',  
                    color_continuous_scale=px.colors.sequential.Plasma, 
                    title='Change in Global Happiness Scores over Available Years')

# Show
choropleth_map.show("notebook")


# Code of 2.3
# List of OECD member countries without abbreviations
oecd_member_countries = [
    'Australia', 'Austria', 'Belgium', 'Canada',
    'Chile', 'Colombia', 'Costa Rica', 'Czechia',
    'Denmark', 'Estonia', 'Finland', 'France',
    'Germany', 'Greece', 'Hungary', 'Iceland',
    'Ireland', 'Israel', 'Italy', 'Japan',
    'South Korea', 'Latvia', 'Lithuania', 'Luxembourg',
    'Mexico', 'Netherlands', 'New Zealand', 'Norway',
    'Poland', 'Portugal', 'Slovakia', 'Slovenia',
    'Spain', 'Sweden', 'Switzerland', 'Turkiye',
    'United Kingdom', 'United States'
]

all_countries = whr['Country name'].unique()  # Get a list of countries
oecd_membership = {country: 1 if country in oecd_member_countries   # Create a dictionary
                   else 0 for country in all_countries}
oecd_df = pd.DataFrame(list(oecd_membership.items()), columns=['Country name', 'OECD'])

# Merge the OECD with the WHR and difference datasets
whr_oecd = whr.merge(oecd_df, on='Country name', how='left')
whr_oecd = whr_oecd.merge(whr_diff, on='Country name', how='left')

# Group by OECD & calculate mean values
numeric_columns = whr_oecd.select_dtypes(include=['float64', 'int64']) # Selecting numeric cols
mean_values = numeric_columns.groupby(whr_oecd['OECD']).mean() 
mean_values_reduced = mean_values.drop(columns=['year'])  # Remove column

# Check if 'year' is in the DataFrame and then drop it if it is
numeric_columns = whr_oecd.select_dtypes(include=['float64', 'int64'])
if 'year' in numeric_columns:
    numeric_columns = numeric_columns.drop(columns='year')

# If 'OECD' is part of the numeric columns, drop it as well since we don't want to standardize it
if 'OECD' in numeric_columns:
    numeric_columns = numeric_columns.drop(columns='OECD')

# Standardize the numeric columns
scaler = preprocessing.StandardScaler()
standardized_values = scaler.fit_transform(numeric_columns)

# Create a DataFrame of the standardized values with original indices and column names
standardized_df = pd.DataFrame(standardized_values, index=numeric_columns.index, columns=numeric_columns.columns)

# Add back the 'OECD' column to the standardized DataFrame
standardized_df['OECD'] = whr_oecd['OECD']

# Melt the DataFrame for visualization
melted_standardized_df = standardized_df.melt(id_vars='OECD', var_name='Metric', value_name='Mean Value')

# Set style of the plot 
sns.set(style="whitegrid")

# Ploting
plt.figure(figsize=(10, 6))
barplot = sns.barplot(x='Metric', y='Mean Value', hue='OECD', data=melted_standardized_df, ci=None)
plt.title('Comparison of Mean Happiness Metrics by OECD Status (Excluding Life Expectancy)')
plt.xticks(rotation=45, ha='right')
plt.ylabel('Mean Value')
plt.xlabel('Metric')

# Show
plt.tight_layout()
plt.show()


# Code of 2.4
from sklearn.model_selection import train_test_split
from sklearn.linear_model import LinearRegression
from sklearn.metrics import mean_squared_error, r2_score
from sklearn.impute import SimpleImputer
from sklearn.ensemble import RandomForestRegressor

# Define predictor variables and target variable
X = whr_na[['Log GDP per capita', 'Social support', 'Healthy life expectancy at birth', 
         'Freedom to make life choices', 'Generosity', 'Perceptions of corruption']]
y = whr_na['Life Ladder']

# Split the dataset into training and testing sets
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, random_state=42)

# Create a Linear Regression model instance
linear_regressor = LinearRegression()
linear_regressor.fit(X_train, y_train)
y_pred = linear_regressor.predict(X_test) # Make predictions
mse_lr = mean_squared_error(y_test, y_pred) # Evaluate the model
r2_lr = r2_score(y_test, y_pred)
# Print the model's coefficients and scores
print(f'For Linear Regression Model:')
print(f'Coefficients: {linear_regressor.coef_}')
print(f'Mean squared error: {mse_lr}')
print(f'Coefficient of determination: {r2_lr} \n')


# Create a Random Forest model instance
random_forest_reg = RandomForestRegressor(n_estimators=100, random_state=42)
random_forest_reg.fit(X_train, y_train)
y_pred_rf = random_forest_reg.predict(X_test) # Make predictions
mse_rf = mean_squared_error(y_test, y_pred_rf) # Evaluate the model
r2_rf = r2_score(y_test, y_pred_rf)
# Print the model's scores
print(f'For Random Forest Model:')
print(f'Mean squared error: {mse_rf}')
print(f'Coefficient of determination: {r2_rf}')


# Actual vs Predicted Happiness Scores using Random Forest
plt.figure(figsize=(8, 5))
plt.scatter(y_test, y_pred_rf, alpha=0.4)
plt.plot([y_test.min(), y_test.max()], [y_test.min(), y_test.max()], 'k--', lw=2)  # Diagonal line
plt.xlabel('Actual')
plt.ylabel('Predicted')
plt.title('Actual vs. Predicted Happiness Scores - Random Forest')
plt.show()


# # Get feature importance
feature_importance = random_forest_reg.feature_importances_
feature_names = X.columns
importance_df = pd.DataFrame({'Feature': feature_names, # Create a DataFrame
                              'Importance': feature_importance})
importance_df = importance_df.sort_values(by='Importance', ascending=False) # Sort in descending order

# Visualizing
plt.figure(figsize=(10, 6))
plt.barh(importance_df['Feature'], importance_df['Importance'])
plt.xlabel('Importance')
plt.title('Feature Importances in Predicting Happiness Score')
plt.gca().invert_yaxis()
plt.show()


# Code of 2.5
numeric_columns2 = whr.select_dtypes(include=['float64', 'int64']) # Selecting only numeric columns
mean_values_counrty = numeric_columns2.groupby(whr['Country name']).mean() # Group by OECD & calculate mean
mean_values_counrty2=mean_values_counrty.dropna().reset_index()
data1=mean_values_counrty2.drop(['year','Country name'], axis=1)
# data1.head()

# Standardization
z_scaler = preprocessing.StandardScaler()
data_z = z_scaler.fit_transform(data1)
data_z = pd.DataFrame(data_z)

# Normalization
minmax_scale = preprocessing.MinMaxScaler().fit(data_z)
data2 = minmax_scale.transform(data_z)
pd.DataFrame(data2).head()

# Calculate the mean distortions for different values of k
K = range(1,11)
meandistortions = []
for k in K:
    kmeans = KMeans(n_clusters=k)
    kmeans.fit(data2)
    meandistortions.append(
        sum(np.min(cdist(data2,kmeans.cluster_centers_,'euclidean'), axis=1)
           ) / data2.shape[0]
    )

# Plot the Elbow Curve
plt.plot(K, meandistortions, 'bx--')
plt.xlabel('k')
plt.show()


# Initialize KMeans with 4 cluster (observed from the Elbow plot)
k_means = KMeans(init='k-means++', n_clusters=4, max_iter=500)
k_means.fit(data2)
label = k_means.fit_predict(data2)

data3 = mean_values_counrty2['Country name']
data4 = data3.unique()

# Check cluster results
data_type = pd.DataFrame(label)
data_type.columns = ['Clusters'] # Name the Column
whr_cluster = pd.merge(mean_values_counrty2, data_type, left_index=True, right_index=True)
# whr_cluster.head()

# Extract countries in cluster 0,1,2,3
cluster_0 = whr_cluster.loc[whr_cluster['Clusters'] == 0, 'Country name'].values
cluster_1 = whr_cluster.loc[whr_cluster['Clusters'] == 1, 'Country name'].values
cluster_2 = whr_cluster.loc[whr_cluster['Clusters'] == 2, 'Country name'].values
cluster_3 = whr_cluster.loc[whr_cluster['Clusters'] == 3, 'Country name'].values

print(f'Cluster 0: \n {cluster_0} \n')
print(f'Cluster 1: \n {cluster_1} \n')
print(f'Cluster 2: \n {cluster_2} \n')
print(f'Cluster 3: \n {cluster_3} \n')


# Boxplots in 3*3 format
fig, axes = plt.subplots(3, 3, figsize=(12, 12))

# List of features
features = ['Life Ladder', 'Log GDP per capita', 'Social support',
            'Healthy life expectancy at birth', 'Freedom to make life choices',
            'Positive affect', 'Generosity', 'Negative affect',
            'Perceptions of corruption']

# Plotting
for i, feature in enumerate(features):
    row = i // 3  # row index
    col = i % 3   # col indext
    sns.boxplot(x='Clusters', y=feature, data=whr_cluster, ax=axes[row, col])
    axes[row, col].set_title(f'Distribution of {feature} by Cluster')

plt.tight_layout()
plt.show()

Name	Variable Description	Type
Country name	Name of the country	Categorical
year	Year in which the data was recorded	Categorical
Life Ladder	Overall happiness score, a measure of well-being	Numeric
Log GDP per capita	Logged Gross Domestic Product per capita, indicating economic status	Numeric
Social support	Degree of social support, based on social networks and individual's feeling of support	Numeric
Healthy life expectancy at birth	Expected number of years of healthy life at birth	Numeric
Freedom to make life choices	Degree of personal freedom in life choices	Numeric
Generosity	Level of generosity, often based on charitable giving assessments	Numeric
Perceptions of corruption	Perceptions of corruption, often based on the level of corruption in public and private sectors	Numeric
Positive affect	Frequency or intensity of experiencing positive emotions or affect	Numeric
Negative affect	Frequency or intensity of experiencing negative emotions or affect	Numeric

Descriptions to the Dataset - World Happiness Report Dataset¶

Research Questions¶

Correlation Between Factors and Happiness Scores:¶

Happiness Trends Over Time in Different Regions:¶

Differences in Happiness Scores Across Economic Blocs:¶

Predictive Analysis of Happiness Scores:¶

Cluster Analysis for Identifying Similar Happiness Traits:¶

Conclusion¶

Code¶